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Cosmological parameter constraints from the CMB power spectra alone suffer several well-known 
degeneracies. These degeneracies can be broken by numerical artefacts and also a variety of physical 
effects that become quantitatively important with high-accuracy data e.g. from the Planck satellite. 
We study degeneracies in models with fiat and non-fiat spatial sections, non-trivial dark energy 
and massive neutrinos, and investigate the importance of various physical degeneracy-breaking ef- 
fects. We test the CAMB power spectrum code for numerical accuracy, and demonstrate that the 
I numerical calculations are accurate enough for degeneracies to be broken mainly by true physical 

(-***) effects (the integrated Sachs- Wolfe effect, CMB lensing and geometrical and other effects through 

recombination) rather than numerical artefacts. We quantify the impact of CMB lensing on the 
power spectra, which inevitably provides degeneracy-breaking information even without using infor- 
mation in the non-Gaussianity. Finally we check the numerical accuracy of sample-based parameter 
constraints using CAMB and CosmoMC. In an appendix we document recent changes to CAMB's 
numerical treatment of massive neutrino perturbations, which are tested along with other recent 

Q\ improvements by our degeneracy exploration results. 
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Observations of the CMB can provide accurate constraints on cosmological models given various fairly weak assump- 
tions. Numerical results for comparison with observations are usually calculated using a linear line-of-sight Boltzmann 
code, with additional modelling of non-linear effects such as CMB lensing and the Sunyaev-Zel'dovich effects. Anal- 
ysis of the WMAP data [T] uses CAMB [2J, a code evolved from an early version of CMBFAST 3 , itself developed 
from the Boltzmann hierarchy code COSMICS @J. As with the other recent codes CLASS [5] and CMBEASY : 6l, CAMB 
and the latest version of CMBFAST aim to calculate the temperature and polarization power spectra at sub-percent 
precision in a matter of seconds, and use of such codes is now routine. Forthcoming data, especially from Planck, will 
greatly increase the precision of available data on small-scales, and it remains important to check numerical robustness 
and understand the physical limitations on what the CMB data can in principle constrain, for example due to the 
well-known geometrical degeneracy. 
\Q Previous work has shown that CAMB is consistent with independent numerical codes and can reliably be used for 

C*") precision cosmology when used with appropriate accuracy settings (and assuming that the physical model including 
recombination history is correct) [7J |5] ■ Older work has also studied the stability of CMBFAST by comparison against 
a full Boltzmann hierarchy code [9]. However these studies have focused on spatially- flat or otherwise restricted 
models. Recent new approximations and coding development aimed towards Planck data analysis have led to faster 
t-H high accuracy calculations, and accuracy settings have been adjusted to allow efficient calculation of power spectra 
at the < 0.1% accuracy level [TPl - fTS] . As the code has developed the approximations and accuracy settings have 
been extensively tested as part of routine development, and we shall not labour the reader by giving an extensive 
presentation here. In this paper we focus on providing a powerful check on numerical accuracy for non-flat (and 
flat) models by exploring parameter degeneracies, where small numerical errors could potentially lead to spurious 
degeneracy breaking. 

This paper has three objectives: (i) testing the numerical accuracy of CAMB and determining the precision that 
is required; (ii) examining to what extent degeneracies are broken by physical processes when precision CMB power 
spectrum data is available; and (iii) demonstrating the extent of degeneracies that can be expected from Planck 
when using CMB power spectrum data alone, independently of parameter priors. We use the numerical code CAMB 1 
exclusively since at the time of writing CLASS [5] does not support non-flat models, and CMBFAST [3] and CMBEASY [13] 
are no longer actively maintained. The July 2011 version of CAMB was originally used for much of this work, and 
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we made various adjustments as indicated by the tests of this paper to maintain accuracy; however for consistency 
of presentation the hnal plots are shown from the October 2011 version after we made various accuracy parameter 
tweaks. In Sec. |VII A| we discuss a further minor modification made following the results of this paper that is now 
implemented in the January 2012 CAMB version; the new version only improves the accuracy of the results shown in 
this paper. 

As is well-known, CMB lensing (for reviews see Refs. [14), IT5] ) can break the geometrical degeneracy since the 
lensing deflections are sourced all along the line of sight, and hence are sensitive to both the geometry and growth 
of structure after recombination |16H18| . Lensing does not bias parameter constraints if it is modelled consistently, 
and the effect on the error bars of the power spectrum estimators is also small [T7] (unless considering low-noise 
B-mode observations [E]). Reference [16] has looked at the effect of lensing on the power spectrum on parameter 
degeneracies for flat and open models, and many authors have also considered the additional information available 
by using the lensing three- and four-point functions (typically by performing lensing reconstruction). Here we focus 
on updating the analysis of lensing on the power spectra to the sub-percent precision era and current understanding 
of the cosmological model, allowing for flat, open, and closed models. Note that the unlensed power spectra are 
not observable directly, so this lensing information is inevitably present in any consistent cosmological CMB power 
spectrum analysis. The effect of lensing on the power spectra has recently been used to constrain dark energy using 
only CMB data from ACT [18] and SPT [20] , but future data constraints will require validation to significantly higher 
precision. The more complicated question of how correctly to use non-Gaussian lensing information in combination 
with the information in the lensed power spectra is left for future work. 

We shall assume the recombination history is accurately calculated, as studied in detail by many authors [2"TH2"T] . 
so that it is not a source of bias or uncertainty. We approximate reionization as being fairly sharp, using CAMB's 
standard parametrization where hydrogen reionization and the first reionization of helium happen together |28j . 

This paper is organized as follows. We start in Sec. [IT] by considering the geometrical degeneracy in non-flat ACDM 
models, considering both numerical and physical effects that determine the extent of the degeneracy. In Sec. Ill we 



then restrict to flat models but allow more general dark energy. In Sec. IV we analyse the approximate degeneracy 
that persists even in a standard flat model with a cosmological constant. In Sees [V] and |VI| we then include massive 
neutrinos, for the case of both flat and non-flat models (Appendix[A]describes recent improvements in CAMB's massive 
neutrino modelling which are tested by these sections). Finally in Sec. |VII we look at parameter constraints expected 
from the Planck satellite using a standard MCMC analysis, demonstrating the extent of the degeneracy expected, 
and quantifying the impact of residual numerical errors on parameter constraints. There, we also describe a new 
interpolation scheme that removes the leading numerical artefact shown in the previous sections; we demonstrate that 
this limits numerical biases in parameter constraints to being < 5% of the random error while maintaining acceptable 
numerical speed. 



II. GEOMETRICAL DEGENERACY WITH A COSMOLOGICAL CONSTANT 



A. Calculating the degenerate models 

In this section we test the numerical accuracy of CAMB using the well-known geometrical degeneracies. We also 
look at physical effects that break the degeneracy including CMB lensing, a geometrical averaging effect through 
recombination, and the late-time integrated-Sachs- Wolfe (ISW) effect. 

We first recap the reason for the geometrical degeneracy. A parameter degeneracy effectively describes our inability 
to distinguish certain cosmological parameter combinations, in this context through using CMB anisotropies alone. 
The primary CMB anisotropies are generated around recombination and what we observe is a projection of conditions 
on the last-scattering surface. If we keep the physical densities in baryons, cold dark matter and the number of 
(masslcss) neutrinos fixed, the pre-recombination physics is unchanged. Since the mapping of physical scales at 
last-scattering to observed angular scales depends only on the angular-diameter distance to last-scattering, there are 
generally degenerate combinations of "late-time" parameters (such as the curvature parameter Q,k and expansion rate 
today Hq) that yield very nearly the same power spectra of primary anisotropies [29-31]. We will consider various 
combinations of parameters that can give nearly identical unlensed CMB power spectra, starting with the degeneracy 
in non-flat ACDM models. 

The first step is to determine which sets of parameters give these nearly identical CMB power spectra. We define a 
fiducial model, and then explore other combinations of parameters that are (nearly) degenerate. We use as fiducial pa- 
rameters the best-fit seven-year WMAP parameters [T] given in Table|lj We run CAMB at boosted high accuracy settings 
to ensure that the fiducial model itself has minimal numerical error: we set the high_accuracy_def ault parameter 
to 'true', and set the additional three accuracy parameters ISampleBoost, lAccuracyBoost and AccuracyBoost to 



3 



Parameter 


Value 




0.733 


n b h 2 


0.0226 


n c h 2 


0.112 


U K 


0.0 




0.0 


h 


0.71 


A s 


2.1364xl0~ 9 


n s 


0.96 


T 


0.088 



TABLE I: Parameter values for our fiducial model. Note that A s is defined at a pivot scale of kg = 0.05Mpc 1 and 

running of the spectral index n s is assumed to be zero. 



values of 2, thus ensuring all numerical calculations are performed with extremely high precision 2 . Henceforth it 
should also be noted that the term 'high accuracy' corresponds to the high_accuracy_def ault parameter having a 
value of 'true' (designed to be appropriate for Planck analysis), whilst 'low accuracy' corresponds to that parameter 
having value 'false' (appropriate for WMAP analysis). 

Once the fiducial model has been specified it is possible to determine which combinations of parameters create 
the geometrical degeneracy by direct comparison of the power spectra returned by CAMB. We want to quantify the 
difference between the spectra in a way that is relevant for observations: for example cosmic variance means that 
much lower accuracy is required on large-scales than on small-scales. What we observe is the CMB temperature T(n) 
along the line of sight n, which can be decomposed into spherical harmonics as 

T(n) =J2 T imYim(n). (1) 

i ni 

The power spectrum for statistically-isotropic fluctuations is defined by 

\TimTii m i) = Su'S mm 'Ci, (2) 
which in a perfect experiment can be estimated by 



Q = 2l 



^yEi^i 2 - ( 3 ) 



From these estimators the probability of the true power spectrum C; (with a flat prior) is given by the likelihood C, 
where for convenience we define x 2 ff = — 21og£. For Gaussian perturbations 



teff 

In 



xle ■ 



In any particular realization of the sky, the estimators C/ have cosmic variance about the true power spectrum, and 
hence scatter in a realization-dependent way. To explore the degeneracies in a realization-independent way we replace 
Ci with the power spectrum in the fiducial model, Cf d . We can then quantify how close a degenerate model is to 



2 ISampleBoost changes the sampling in I over which the C; are interpolated. lAccuracyBoost changes the number of multipoles maintained 
when integrating the Boltzmann hierarchies. AccuracyBoost changes wavenumber sampling and integration step sizes (and various 
other parameters). The high_accuracy_def ault switch increases the density and range of wavenumber samples, increases the number of 
multipoles retained in the Boltzmann hierarchies that are evolved, switches from the tight coupling approximation slightly earlier, uses a 
larger range of unlensed template C; when calculating the lensed C;, and, prior to January 2012, increases the density of £ samples that 
are interpolated (see Sec. |VII At . Unlike blindly increasing the accuracy parameters, which changes many parameters simultaneously 
by the same amount, the changes adopted by the high_accuracy_def ault switch give the minimal internal parameter tweaks required 
to give the target accuracy < 0.1% at I > 500 and hence do not dramatically increase the running time. For further discussion see 
Refs. 0|8], though changes have been made since those papers. 
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our fiducial model by using the effective chi-squared value \cfi- Note that we have normalized this equation such that 
any model that is exactly equal to the fiducial model will give a value of \cb = 0- Throughout, we use a value of 
^max = 2000 with Z min = 2 or l min = 1 00 as indicated below (depending on whether we are focusing on high-Z numerical 
errors or a more realistic analysis). The cut at i max = 2000 is somewhat arbitrary, but reflects that lower precision is 
required on small-scales due to rapidly growing uncertainties from beams, point sources, secondary anisotropies, and 
other foregrounds. In practice the numerical performance of CAMB above / max = 2000 is not dramatically worse, and 
should be sufficient if data on smaller scales can be used. 

By minimizing Xcfii we can nno - se ^ s °f parameters that give power spectra that are very close to the fiducial model. 
Example unlensed power spectra are shown in non-flat ACDM models in Fig.[l] In this case, the geometric degeneracy 
is within the two-dimensional space of Q\ and h. 




10 100 1000 10 100 1000 

L (multipole) / (multipole) 

FIG. 1: CMB power spectrum obtained using CAMB for nearly degenerate geometries in non-flat ACDM models 
with no lcnsing (left) and the fractional differences from the fiducial-model spectrum (right). Both VLbh 2 and fl c h 2 
were fixed to their fiducial values in all cases to preserve the pre-recombination physics. Low accuracy and values of 
1 for ISampleBoost, AccuracyBoost and lAccuracyBoost were used for the calculations. 



B. Degeneracy breaking effects 

1. Numerical accuracy 

As shown in Fig. [T] the power spectra are not entirely degenerate. On large-scales this is due to the late-time ISW 
effect; on small scales the power spectra are much closer. However the small-scale spectra will not be quite identical 
due to both numerical artefacts and small physical effects. For the purpose of having reliable parameter constraints 
we would of course like the numerical effects to be small, so the degeneracy is broken only by the physical effects. 
Here, we quantify the difference between these degenerate power spectra as a check that the numerical accuracy of 
CAMB is sufficient for forthcoming data, and investigate in more detail degeneracy breaking effects such as lensing, 
geometrical effects through last scattering, and also the late-time ISW effect. 

We expect any numerical errors in the power spectrum as a result of CAMB's calculations to be comparable to, 
if not less than, the quoted values: 0.3% for low accuracy calculations (high_accuracy_def ault=F) and 0.1% (at 
I > 500) for high accuracy calculations (high_accuracy_def ault=T). Figure]!] (right-hand panel) shows the fractional 
differences between the power spectra in four degenerate geometries and the fiducial model. The spectra are closely 
similar at high I where the ISW contributions to the spectra all become small. 

Figure [2] shows differences between the high-^ power spectra for a selection of nearly-degenerate models. We can see 
a larger periodic oscillation in the difference plots that roughly corresponds to the location of the peaks and troughs in 
the power spectra. These are dependent on the degenerate model used (and increase with |S1_r-|), which indicates they 
are due to physical effects (see discussion in the next subsection). However, we also see much smaller more random 
oscillations within this larger periodicity. These are due to numerical errors and are reduced slightly by increasing 
from low accuracy calculations to high accuracy calculations; we shall quantify whether this improvement is good 
enough in more detail later. Boosting the accuracy removes most of the residual numerical error, and we are left with 
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smoothly varying differences due to physical effects. 

We can nearly isolate the numerical errors in the computation of a given model by subtracting a very accurate 
computation (high accuracy and all parameters boosted to 2) of the same model. Our results are plotted in Fig. [3] 
We see that the numerical error for low accuracy calculations is of the order of 0.2% (rather better than the quoted 
accuracy, indicating that the calculation is slower than it needs to be) . There is a small increase in numerical precision 
when high accuracy calculations are used, decreasing the error to within the quoted value of 0.1% at I > 500 where 
cosmic variance becomes small. 
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High Accuracy, Parameters = 1 




High Accuracy, Parameters = 2. 
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FIG. 2: Difference between the unlensed power spectra of a range of non-flat degenerate models and the fiducial 
model. Models are computed at low accuracy with default accuracy parameters (left), high accuracy with default 
accuracy parameters (middle) and high accuracy with parameters boosted to 2 (right). The right-hand figure shows 
the physical differences in the spectra, with very little residual numerical error. 
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FIG. 3: Numerical errors in the computation of a range of non-flat degenerate models. For a given model, the 
numerical errors are estimated by subtracting the spectrum from one calculated at high accuracy with accuracy 
parameters boosted to 2. Errors are plotted for a low (left) and high (right) accuracy calculation. In all cases, the 

other accuracy parameters are at their default values. 
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If numerical precision is the only objective, increasing the accuracy settings will give accurate smooth results. 
However, such settings can increase the running time required very significantly, so rather than simply increasing 
these parameters it is best to optimize them such that a balance between speed and precision is found (see Ref. [7] for 
previous work on this). It is easy to identify the main cause of the high-frequency numerical wiggles: this is just due 
to small interpolation errors from the rather sparse I sampling used when calculating the power spectra (Al ~ 50). 
This effect is suggested by the clear correlations between models of the numerical errors isolated in Fig. [3] Figure [4] 
shows that increasing ISampleBoost to 2 can effectively remove these wiggles, as expected. However using twice 
the density of I samples takes nearly twice as long (in non-flat models) to calculate spectra. Amongst other things, 
setting high_accuracy_def ault lowers the I sampling at high I to Al = 42 from = 50, which significantly lowers 
the residuals; however this still leaves some interpolation artefacts. A further boost of ISampleBoost ~ 1.2—1.5 is 
enough to remove most of these if required (so that Al — > AZ/lSampleBoost) while maintaining good efficiency. A 
better alternative to increasing blindly the I sampling is to improve the interpolation method as discussed later in 
Sec. VII A On large scales the numerical errors are larger, ~ 0.3%, but still well below the physical differences from 
different late-time ISW contributions between nearly-degenerate models, and always smaller than cosmic variance. 




FIG. 4: Fractional differences between nearly degenerate geometries and the fiducial model, comparing boosts of 2 
in all the accuracy parameters (solid red) to boosts of 2 in only ISampleBoost (dashed blue). In all cases, high 

accuracy is used for the calculations. 




FIG. 5: Minimum x^g for a range of degenerate geometries close to the fiducial model. Solid (black) is high accuracy 
and each of ISampleBoost, lAccuracyBoost and AccuracyBoost boosted to 2; dashed (green) is high accuracy with 
default accuracy parameter values; dot-dashed (blue) is high accuracy with only ISampleBoost boosted to 2; and 
triple-dot-dashed (red) is low accuracy with no accuracy boosts. In all cases we use only I > 100 in the calculation 

01 Xcff ■ 



We now test to what extent numerical errors break the geometric degeneracy compared to physical effects by 
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comparing, for a range of values of J7a, the minimum values of Xcs (with the fiducial model for C/ and including only 
multipoles I > 100) with respect to h. Figure [5] shows the minimum \cS as a function of close to the fiducial 
value; the models are calculated for a range of accuracy settings. The solid black line indicates the closest we can get 
to our fiducial model in terms of numerical accuracy, and the shape of the curve reflects almost entirely the physical 
effects that weakly break the degeneracy; these are quantitatively more important the further away the degenerate 
geometry is from the fiducial model. 

Numerical errors change this result, shifting the curve to a different offset and shape. The low accuracy settings 
are significantly biased (by around 0.5cr) but use of default high accuracy settings recovers a much more accurate 
likelihood curve. However the likelihood curve is still slightly shifted compared to the accurate calculation. If we 
increase the value of ISampleBoost to 2, effectively removing all the small-scale interpolation wiggles, we find the 
curve returns to nearly the correct shape, with only a slight offset due to residual numerical effects. Offsets do not 
affect inferences about parameter constraints and are relatively harmless. However shifts in the curves give rise to 
errors in the posterior mean and maximum likelihood point, and give rise to a parameter biases. We assess this in 
more detail below (Fig. [9]) in the more observationally-relevant case that both noise, \ow-l modes and lensing are 
included. 



2. Physical effects and Planck data 

Figure |2] shows that the geometric degeneracy for non-flat models is broken at the 10 -3 level (for flxh 2 ~ 0.02) 
by physical (geometrical) effects even in the absence of lensing. The effect of curvature on the evolution until 
recombination is very small, 0(1O -4 ), so the physical anisotropy sources are essentially the same in the various 
degenerate models. However there is a larger effect due to the finite thickness of the last-scattering surface (c.f. 
Ref. [32 ). The change in the (comoving) angular-diameter distance d\ over a fixed thickness of last-scattering 
centred on a given (1a{z*) (fixed by the peak of the visibility function) in curved models differs from the flat case; it 
is smaller in closed models and larger in open. In detail, in a closed model <Ia — sin(v / ^x)/V r ^ at radial distance 
X, so that dd,A/dx = cos{y/K\) = (1 — Kd A ) x / 2 for \[K\ < it/2. In an open model, the equivalent result is 
ddA/dx = (1+ \K\d 2 A ) 1 / 2 . If we look at an angular scale corresponding to an acoustic peak — probing a perturbation 
scale that is a maximum or minimum at the peak of the visibility — then we will also generally see slightly larger scales 
at earlier times near the start of recombination, which are not yet at an extremum. However, in a closed universe the 
change in transverse comoving scale going through last-scattering is smaller: the perturbations probed earlier would be 
closer to an extremum. In an open universe there is therefore a slight suppression of the power in the acoustic peaks, 
and in closed universe a slight boost, as shown in Fig. [2] Assuming the thickness of recombination is a significant 
fraction of 77*, the conformal time at the peak of the visibility function, this effect is 0(r}*dAK) ~ 10 -3 -10~ 2 for the 
parameter range of interest. This is large compared to other effects of curvature on the pre-recombination dynamics 
because we are viewing last scattering from a large distance ( dA 3> t]*), and the effect of curvature on the fractional 
change in dA through a fixed thickness increases with radial distance. 

In principle perfect small-scale unlensed data could therefore determine the curvature, even without additional 
information. However in the immediate future, the best data available on intermediate scales will come from the 
Planck satellite. To assess more carefully the accuracy required for Planck we need to look in more detail at what is 
actually observed, and the relevant errors. For example the CMB is inevitably lensed, so we only actually observe the 
lensed power spectra, and it is the accuracy of parameter constraints including the lensing effects that really matter. 
We also introduce a very simple noise in the form relevant for Planck. 

In terms of the power spectrum, lensing serves to smooth out the peaks corresponding to the CMB anisotropics [16] . 
For example, random gravitational lenses between the observer and the surface of last scattering cause the temperature 
anisotropies that we see effectively to be 'smeared' out, causing us to loose definition in these regions and damping 
the sharpness of the features. This gives a smoothing of the power spectrum at the peak locations. This causes a 
breaking of the geometric degeneracy, as the amount of lensing is dependent on the geometry and growth of structure 
between us and last scattering [TMISl [20] . 

When we include lensing the differences between the previously nearly-degenerate geometries increases significantly, 
as shown in Fig. [6} Whilst the periodicity we observed earlier is still present, the amplitude of these peak differences 
increases nearly ten-fold on small scales as the different geometries cause different amounts of lensing, and hence 
different amounts of smoothing in the power spectra. As we will see later this is also apparent quantitatively, where 
we see a large increase in the minimum effective chi-squared values for each degenerate model in the presence of 
lensing, allowing us to place a much better constraint on tt\ (or fix)- 

Including noise in the procedure, on the other hand, has the opposite effect to lensing, in that it strengthens the 
degeneracy by decreasing the sensitivity to small changes in the power spectra at high I. As a simple test we include 




FIG. 6: Fractional change in the lensed CMB power spectra at different nearly-degenerate non-flat models. Note 
that the fractional differences between power spectra are much larger with lensing included, so here we plot a 
significantly narrower range of nearly-degenerate models than in the unlensed case of Fig. [2] Here results are 
calculated at high accuracy with ISampleBoost = 1.17 to remove most of the small interpolation wiggles. 



noise on the temperature, Nt, and polarization, Ng, spectra with 

N T = ^ = N e l « +1 V° 2 , (5) 

where No = 0.5 x 1CP 4 fiK 2 is the white-noise level (corresponding to 24//Karcmin) and a = 1.7 x 10~ 3 rad for a 
beam of 7arcmin full- width at half-maximum. 

Our choice of noise level is at the lower end of what is likely to be achievable by combining multiple sky scans 
and frequencies from Planck, and is therefore conservative (in that if numerical accuracy is sufficient in our tests, it 
will almost certainly be sufficient in reality). We then simply add this noise to the polarization and temperature C;s, 
Cf E and Cj T respectively, returned by CAMB. Since the polarization signal is much smaller than the temperature it 
has significant noise, but nonetheless can have some degeneracy breaking power (as well as constraining the optical 
depth). For current purposes we can neglect any i3-mode polarization, and the xis including polarization (and noise) 
is then given by 



Xoff 



= £(21 + 1) ('&(C,Cr 1 )+ln 




- 2 



where Ci = 



(jTE (jEE 



(6) 



Again note that this has been normalized so that any degenerate model exactly equal to the fiducial model Ci = Cf d 
will return a value of xls = 0- 

Once noise is included, constraints become significantly weaker, but this is compensated by lensing degeneracy 
breaking. Maximum likelihood curves are shown in Fig. [7J where we now also vary other cosmological parameters, 
as small changes in these may be able to compensate partly degeneracy-breaking due to lensing and other physical 
effects: in addition to h we minimize over fl^h 2 ,tt c h 2 , A s ,n s and r for each value of Ha- As expected, even in the 
presence of noise the effect of lensing is enough to break the degeneracy and constrain Q\ well. When we use a 
six-parameter fit the breaking of the degeneracy is lessened somewhat, though not substantially, indicating that the 
lensing effect cannot easily be mimicked by other changes of parameters (at least if dark energy is assumed to be a 
cosmological constant). 

We also see from Fig. [7] that the late-time ISW signal in the low-Z anisotropics can break the degeneracy significantly, 
though the lensing information makes the constraint about a factor of two tighter. This is shown further in Fig. [8] for 
a full MCMC parameter analysis using CosmoMC [33] from the same likelihood but with four times larger noise (as 
might be nearer to reality if only one frequency is used). 

We finally assess the impact of the numerical errors we found previously, in particular the small shifts in the inferred 
value of f^A when different values of ISampleBoost were used. Figure [9] shows the likelihood curves for various values of 
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FIG. 7: Minimum effective chi-squared values as a function of We consider unlensed spectra with noise, using 
I > 2 (red triple-dot-dashed) or only I > 100 (black solid), and lensed spectra with noise using I > 2. With lensing 
included, we consider both a simple minimisation with respect to h (green dashed) and a standard six-parameter fit 
(Multifit) minimising with respect to h, f^ 2 , Q c h 2 , the spectral index, n s , the optical depth, t, and A s e~ 2r (blue 

dot-dashed). 



ISampleBoost with both noise and lensing included. High accuracy with default settings is slightly shifted compared 
to a more accurate calculation, though the shift is not large compared to the overall degeneracy width. For a more 
accurate calculation ISampleBoost > 1.2 could be used. 



III. GEOMETRICAL DEGENERACY IN A FLAT UNIVERSE WITH A MORE GENERAL DARK 

ENERGY 



In this section we investigate the geometrical degeneracy in a flat universe when there is dark energy parameterized 
by a constant equation of state parameter w = P/p (we restrict to w < —1/3). We take a fiducial value w = — 1 
corresponding to a cosmological constant. Figure [lO] shows four different nearly degenerate models obtained by 
varying w and h. Very different values for w can still produce nearly identical CMB anisotropics: the degeneracy is 
only significantly broken at low I by the late-time ISW effect. The flat-model dark energy degeneracy is more exact 
in the unlensed CMB than the non-flat case, and hence provides a stringent test of numerical accuracy of CAMB as 
the geometry is varied (though a much less sensitive test of errors in the pre- recombination physics) . We focus on the 
high-Z regime as an accuracy test, where differences in the spectrum are well below a percent. 



A. Degeneracy breaking effects 



1. Numerical Accuracy 



Figure 11 shows the change in the unlensed power spectrum for different degenerate models at various accuracy 
settings. As expected we see some small numerical wiggles at standard settings which subsequently disappear at 
boosted accuracy. As expected there is hardly any breaking of the degeneracy due to physical effects at I > 100. The 
degeneracy is much more accurate than in the non-flat case (with w = — 1) because the physical effect of the dark 
energy is negligible in the early universe, and the distance to last-scattering and geometry is now fixed, so there is no 
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FIG. 8: Idealized forecast parameter constraints from Planck power spectra only, varying A s , n s , ttbh 2 , £l c h 2 , r, 9, 
and VLk with flat priors, and Hq, fi m and J7a being derived parameters. Points show samples from the expected 
posterior if the unlensed power spectra were observed, the black contours the better constraint obtainable in reality 
accounting for power spectrum lensing. Adding lensing reconstruction information could further shrink the extent of 

the degeneracy by a factor of roughly two. 
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FIG. 9: Minimum effective chi-squared values for non-flat geometries, in the presence of both lensing and noise, 
using high accuracy and varying the value of ISampleBoost. Note that here we use alH > 2 to make this figure 

comparable to Fig. [7] 

change in the geometrical effect from averaging through last-scattering. 

At low accuracy, the numerical errors are, on the whole, less than the quoted accuracy of 0.3%. There is a small 
'drift' at low accuracy; however, even including this drift, the errors are still generally significantly below 0.3% at high 
I. The small 'drift' disappears with high accuracy calculations. 

We can quantify how much these numerical errors affect the w-h degeneracy by looking at the minimum effective 
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FIG. 10: Geometric degeneracy in flat models with w ^ — 1. Unlensed CMB power spectra for several different 
nearly degenerate flat models are plotted (left), along with their fractional differences with respect to the fiducial 
model (right). The fiducial model is the same as used in Sec. ITT] (and so has w = —1). 
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FIG. 11: Fractional differences between unlensed power spectra and the fiducial model for flat dark energy models. 
The nearly-degenerate models are computed at low accuracy with all accuracy parameters set to 1 (left), and at high 
accuracy with all accuracy parameters set to 1 (middle) and 2 (right). With sufficient numerical accuracy, the high-Z 
unlensed power spectra are essentially identical for different flat dark energy models unless w is very close to —1/3. 



chi-squared values obtained in this space for various settings of the accuracy parameters. Figure [12] plots the minimum 
effective chi-squared in the w-h space for —2.0 < w < —0.35 at both standard low and high accuracy settings as well 
as boosted accuracy. For standard accuracy there is an offset, in part because of the small numerical 'wiggles' shown 
in Fig. |11| However an offset is usually harmless as it does not affect relative parameter constraints, as demonstrated 
in the figure by renormalizing the curves to have global minimum values equal to zero. At low accuracy the curves 
are not very smooth also due to numerical effects, though only at a low level that goes away when high accuracy is 
used. The small residual shift in shape compared to boosted accuracy settings can be mostly eliminated as in the 
previous section by a slight increase of ISampleBoost to 1.1—1.2. 
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FIG. 12: Left: Minimum effective chi-squared values in the w-h space for flat dark energy models, using all I > 2. 

The models are computed at low accuracy with accuracy parameters set to 1 (red, triple-dot-dashed), and at 
high-accuracy with accuracy parameters set to 1 (blue, dot-dashed) and 2 (black, solid). In all cases, the fiducial 
model has w = — 1. Right: Rcnormalized chi-squared curves so that the global minimum value is equal to zero. 
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FIG. 13: Minimum effective chi-squared values as a function of w in flat dark energy models. We consider unlensed 
spectra with Planck-like noise, using I > 2 (red triple-dot-dashed) or only I > 100 (black solid), and lensed spectra 
with noise using I > 2. With lensing included, we consider both a simple minimisation with respect to h (green 
dashed) and a standard six-parameter fit (Multifit) minimising with respect to h, flbh 2 , £l c h 2 , the spectral index, n s , 

the optical depth, r, and A s e~ 2T (blue dot-dashed). 



2. Effect of noise and lensing 



We now use temperature and -E-mode polarization, as in Eq. [6j to explore how lensing, noise and performing a 
six-parameter fit affect the geometric degeneracy in flat models with general dark energy (w =/= —1); see Fig. 13 There 
is weak degeneracy-breaking at low / due to the late-time ISW effect, as well as additional breaking from the lensing. 
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In the absence of lensing and low I information we see that the degeneracy is almost exact. Finally we can see that 
fitting Qbh 2 , £l c h 2 , A s e~ 2r , n s and r as well as h, when using lensing, gives us a slightly better minimization, but the 
effect is small so that the allowed changes in other parameters cannot effectively mimic the effect of lensing on the 
power spectrum. 



IV. ACOUSTIC-SCALE DEGENERACY IN A FLAT ACDM UNIVERSE 




FIG. 14: Top: Power spectra for four different flat ACDM models that are all nearly degenerate with the fiducial 
model. Bottom: Fractional difference between these degenerate models and the fiducial model. 



In this section we investigate one further CMB parameter degeneracy, mostly involving il m and H , in flat ACDM 
models. This is not a purely geometric degeneracy — it involves multiple parameter variables — but is an important 
source of uncertainty in individual parameter constraints from the CMB alone [2*0 OH]. Depending on the data 
used, ncarly-dcgcncratc models have fl m h 3 approximately constant. The exact direction of the degeneracy varies 
slightly with Z max and noise under consideration, but degenerate models are such that the observed angular scale 
of the acoustic peaks is nearly constant. In a flat ACDM model this can only be achieved by changing both the 
sound horizon at recombination and the angular diameter distance. However changing the sound horizon also involves 
changes in matter densities, which have other effects on the observed power spectrum: the degeneracy only exists to 
the extent that changes in other parameters can compensate these effects within the limits of cosmic variance and 
observational noise. 

For all the results given in this section we maximize the likelihood over f2&/i 2 ,fl c h 2 , A s e~ 2r ,n s and r. We use the 
same fiducial model as before (Table [jj which has ft m = 0.267. For each value of tt m we minimize xls with respect 
to the other parameters. 

Figure 14 shows the power spectra of four nearly degenerate models, spanning a wide range of £l m . The largest 
fractional differences are at low I where the late-time ISW effect breaks the degeneracy due to the different late- 
time expansion histories. However, there are also significant differences between the spectra on small scales so with 
precision data these models can be relatively easily distinguished. However the approximate degeneracy will remain 
over a smaller range of parameter values. 

Figure [15] shows the differences between the power spectra and the fiducial model for nearly-degenerate models over 
a small range of il m about the fiducial value. The differences seen in high accuracy calculations with boosted accuracy 
parameters are due to physical effects. Unlike in previous sections, the differences vary in a complicated way across 
the whole spectrum, stemming from differences in the relative acoustic peak amplitudes and tilt as all the parameters 
are changed (not quite managing to cancel completely the effect of the changing matter densities at recombination). 
There are of course also numerical artefacts which are clearly visible at low accuracy, but largely disappear when 
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FIG. 15: Fractional differences between the CMB power spectra of nearly-degenerate models and the fiducial model. 
The former are calculated using low accuracy and default accuracy parameters (left), high accuracy and default 
parameter accuracy parameters (middle) and high accuracy and boosted accuracy parameters (right). The 
right-hand panel leaves mainly the physical differences between the models. 
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FIG. 16: Numerical errors in the computation of flat ACDM models that are nearly degenerate with the fiducial 
model. For a given model, the numerical errors are estimated by subtracting the spectrum from one calculated at 
high accuracy with accuracy parameters boosted to 2. Errors are plotted for a low (left) and high (right) accuracy 
calculation. In all cases, the other accuracy parameters are at their default values. 



we use boosted high accuracy settings. Following our earlier treatment of geometric degeneracies, we can isolate the 
numerical effects by considering the difference between the model spectrum and one calculated at high accuracy with 
accuracy parameters boosted to 2. Figure [16] shows that the numerical errors at intermediate and high I are below 
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FIG. 17: Left: Minimum effective cm-squared values as a function of J7 m in fiat ACDM models close to the fiducial 
model. The models are computed at low accuracy with default accuracy parameters (red, dot-dashed), at high 
accuracy with default accuracy parameters (blue, dot-dashed) and at high accuracy with accuracy parameters 
boosted to 2 (black, solid). In each case, only temperature data is used and there is no instrument noise. Right: 
Impact of including polarization data with no noise (black, solid), with Planck-like noise (red, triple-dot-dashed) and 
with lensing and noise (blue, dot-dashed). The models are computed at high accuracy with ISampleBoost increased 

to 1.17. 




FIG. 18: Effect of noise in locating nearly-degenerate flat ACDM models by chi-squared minimisation. Fractional 
differences from the fiducial model are plotted for four nearly degenerate models; the latter are computed at high 
accuracy and with ISampleBoost = 1.17 to remove most of the residual interpolation wiggles. Compared to the 
models in Fig. 15 larger differences in the spectra are allowed at high I where noise dominates. 



the quoted accuracies of 0.3% (low accuracy) and 0.1% (high accuracy). 

Although these numerical inaccuracies are within the quoted ranges, Fig. [17] shows that the numerical errors can 
still give a very small shift in the likelihood. As with the similar results for the non-flat fl\-h degeneracy, increasing 
ISampleBoost slightly removes some of this effect. Including noise-free polarization data helps to break the degeneracy 
considerably but lensing has little effect. For Planck-like noise, the degeneracy is broadened slightly from the noise- 
free temperature-only case: noise allows somewhat larger differences in the power spectra at high I where the noise 



becomes significant, as shown in Fig. 18 



Although the degeneracy is quite tight, it does not extend over a large range of parameter values; as such including 
lensing barely has any effect on the chi-squared values returned for each model. In Fig. [19] we show how the best-fit 
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FIG. 19: Best-fit values of Ho as a function of fl m without lensing or noise included (black, solid), with only lensing 
(green, dashed), with only noise (red, triple-dot-dashed), and with lensing and noise (blue, dot-dashed). Optimized 
accuracy settings were used for all calculations. The degeneracy is close to ft m h a ~ const., with 3 < a < 3.1 having 

slightly different values in the various cases. 

value of Hq varies as a function of il m . We can see that lensing has a negligible effect on the degeneracy direction 
because the parameters are so strongly constrained without it anyway. However, noise does slightly change the 
direction; this is expected since the relative importance of different physical effects changes as a function of I. 

V. GEOMETRICAL DEGENERACY IN A FLAT UNIVERSE WITH MASSIVE NEUTRINOS 

In this section, we consider an approximate geometrical degeneracy that arises in models with massive neutrinos 
that are light enough (masses well below 1 eV) to still be relativistic at recombination. We give a brief description 
of CAMB's numerical calculation with massive neutrinos in Appendix. [XJ including recent changes giving improved 
performance that are tested for consistency by the analysis here. 

We consider a flat fiducial model with the parameters given in Table [I] but with the physical density in massive 
neutrinos il^h 2 = 6.0 x 10~ 4 (and hence a lower J7a to preserve flatness). We calculate the fiducial model at high 
accuracy and with all three accuracy parameters (lSampleBoost, lAccuracyBoost and AccuracyBoost) set to 2 so 
the numerical errors are very small. Throughout this section, we assume three neutrino mass eigenstates with squared 
mass differences given by neutrino oscillation results in the 'normal hierarchy'. In particular, we use the central values 
from the 2006 update of the global fit in Ref. [35] (also adopted in Rcf. (36j): 



with mi the lightest neutrino mass, and the heaviest. While these central values are only consistent with more 
recent fits, e.g. Ref. |37j , at around the 2a (of the marginal errors) this should not impact our conclusions on numerical 
robustness and physical breaking of degeneracies. The energy density of neutrinos in our fiducial model corresponds 
to the minimal-mass (mi = eV) normal hierarchy. 

We consider the two-parameter approximate degeneracy between Q, v h 2 and h in flat models. All other parameters 
are fixed to their fiducial values throughout this analysis (including fl c h 2 ; only changes to preserve flatness). 
Physically this degeneracy arises because, on small scales, the only significant change to the power spectrum is through 
a change in the angular diameter distance to last scattering, cLa{z*). For the sub-eV neutrino masses considered here, 
the neutrinos are relativistic at recombination, so the difference they make on the dynamics of the pre-recombination 
universe compared to massless neutrinos is very small. However, at late times at least two of the mass eigenstates 
become non-relativistic increasing the energy density relative to a massless model with otherwise the same parameters. 
This decreases (Ia{z*) and causes the CMB acoustic peaks to shift to larger angular scales, since the sound horizon 
at recombination has only a weak dependence on (light) neutrino masses. The value of (1a(z*)> and hence the peak 
locations, can be restored by reducing the Hubble constant (and hence decreasing the dark energy density which is a 
derived parameter in our analysis of flat models) . 




\ = 7.9 x 10~ 5 eV 2 , 
\ = 2.2 x 10~ 3 eV 2 , 



(7) 
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FIG. 20: Fractional differences between four degenerate models and the fiducial model. The degenerate models are 
computed at low accuracy with no boost in accuracy parameters (left), high accuracy with no other accuracy boost 
(middle) and at high accuracy with the accuracy parameters boosted to 2 (right). The right-hand panel reveals 
small physical degeneracy-breaking effects between the models. 



We construct several nearly degenerate models by locating parameters in the Q, v h 2 -h plane all with the same ratio 
of the sound horizon at last-scattering to cLa{z*) as in the fiducial model. Examples of the fractional differences 



between the unlensed power spectra of these models and the fiducial model are shown in Fig. 20 The most massive 



degenerate model has summed masses of 0.38 eV corresponding to almost degenerate masses mi sa rri2 = 0.122 eV 
and m,3 = 0.131 eV. In the right-hand panel of Fig. [20j the degenerate models are calculated at high accuracy and 
with all accuracy parameters boosted to 2. The degeneracy is nearly exact on small scales with physical differences 
in the power spectra below 0.2% for 100 < I < 2000; the difference would be only marginally detectable (assuming a 
known template shape for AC; and all other parameters fixed) for the our largest mass case with perfect temperature 
data. 



A. Degeneracy breaking effects 

1. Numerical accuracy 



We see the same general trends in Fig. 20 to variation of the accuracy settings as in the models discussed earlier. 
At low accuracy, there is high frequency numerical noise and slowly-varying numerical drifts in the spectra. However, 
these numerical effects are comfortably smaller than the quoted accuracy of 0.3%. With high_accuracy_def ault set 
true, most of the drifts are removed and the amplitude of the high-frequency noise is significantly reduced though not 
eliminated entirely. The high-frequency noise can be further reduced by boosting ISampleBoost or, better, with the 
improved interpolation method introduced in Sec. |VII A| 

We now compute XeS fr° m Eq. (kD) up to i max = 2000 along the degeneracy line in the £l v h?-h plane. Results are 



plotted in Fig. 21 for a range of accuracy settings. The magenta curve is computed with the same high accuracy as the 
fiducial model used as input to xts an d therefore represents our best estimate of the truth. In this case, the minimum 
XcH ls achieved at the fiducial value fl u h 2 — 0.6 x 10~ 3 , and the slow variation of xls ar i ses from the degeneracy not 
being physically exact (see later). Comparing to the low accuracy calculations and default parameter settings, we see 
that in this case the shape of the likelihood is significantly altered with the minimum of Q, v h 2 displaced from the true 
value by ~ la along the degeneracy direction, and the width of the distribution artificially tightened. These errors 
likely arise from the low frequency numerical errors visible in Fig. |20| 



We further see in Fig. 21 that default high accuracy setting get close to the 'true' likelihood shape (with a slight 



offset). The numerical noise that remains in this case (see Fig. 20 1 has only a small effect on the likelihood, and is 
removed by a further boost in accuracy. 
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FIG. 21: Effective mean \ 2 about the fiducial model for a range of degenerate models and accuracy settings: low 

accuracy with default values (red) or boosted values of 2 (blue) for ISampleBoost, lAccuracyBoost and 
AccuracyBoost; and high accuracy with default (green) and boosted (to the value 2; magenta) parameters. In all 
cases, the fiducial model in is computed at high accuracy and all parameters boosted to 2. Only unlensed 

temperature spectra are included and there is no noise. 



2. Physical effects 



The degeneracy is not exact as illustrated in the right-hand panel of Fig. [20] On the large scales relevant for the late- 
ISW effect, the neutrinos are able to cluster and their additional contribution to the expansion rate at late time does 
not impede the growth of structure. However, simultaneously reducing the dark energy density to keep cLa{z*) fixed 
lessens its tendency to halt structure formation and so there is less decay of the gravitational potential and a smaller 
late-ISW effect. On intermediate scales near the first peak, the early-ISW effect is significant. This arises from the 
decay of the gravitational potential as the equation of state changes during the matter-radiation transition. Compared 
to the massless case, neutrinos of mass m„ that are still relativistic at temperature T v (z) have their pressure reduced 
by a fractional amount 0(m v /ksT^) and the energy density increased by a similar fraction. This leads to additional 
decay of the gravitational potential shortly after recombination (i.e. when the neutrinos are no longer ultra-relativistic 
but still make a significant contribution to the energy density) and an increase in the early-ISW effect. Further effects 
of the enhanced energy density are small reductions in the physical sound horizon and damping length but the former 
is compensated by adjusting cLa{z*)- Moreover, the dynamics of the small-scale neutrino perturbations themselves 
acquire 0(1 — w^) corrections, where the typical neutrino thermal speed v„(z) = 1 — ©(m^/fc^T^) 2 changes the free- 
streaming scale from the (light) horizon which feeds back into the evolution of other perturbations. The expected 
fractional size of these effects on the CMB power spectrum is ~ 0(rn v / 'ksT^) 2 p v / ptot, where p v is the neutrino energy 
density and p to t the total energy density. In our most massive model, this suggests effects at the O(10~ 3 ) level. The 
right-hand panel of Fig.[20]shows the relative physical differences of our degenerate model power spectra to the fiducial 
model. The expected changes in the late- and early-ISW effects can clearly be seen. In addition, on smaller scales we 
see effects of 0(1O~ 3 ) for our most massive model, consistent with the estimate above. 

Note that our analysis here does not include the effect of weak gravitational lensing of the CMB. If this effect were 
included, the degeneracy would be further broken since massive neutrinos reduce the power spectrum of the lensing 
deflections on small scales compared to a model with massless neutrinos. 
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FIG. 22: Fractional differences between nine degenerate models (three in each row) and the fiducial model exploring 
the J7 y /i 2 -J7if-/i degeneracy. The degenerate models are computed at low accuracy with no boost in accuracy 
parameters (left column), high accuracy with no other accuracy boost (middle column) and at high accuracy with 
the accuracy parameters boosted to 2 (right column). The right-hand panels reveal small physical 
degeneracy-breaking effects between the models. Moving from the top row to the bottom at fixed line style (colour), 
the neutrino mass is increasing at fixed h (and 9). Within a panel, the mass is constant. 



VI. THREE-PARAMETER GEOMETRIC DEGENERACY WITH MASSIVE NEUTRINOS AND 

CURVATURE 

We now extend the parameter space of the previous section to include spatial curvature. We expect a three- 
parameter geometric degeneracy involving h, Q, v h 2 and fix at fixed 9 since neither curvature nor sub-eV neutrino 
masses are dynamically significant before last-scattering. We retain the fiducial model used in the previous section, 
and continue to probe only sub-eV neutrino masses. We do not include lensing. 

We find degenerate models by fixing the angular size 9 of the sound horizon to its value in the fiducial model. The 
likelihood is not exactly constant over the two-dimensional degeneracy surface due to the small physical degeneracy- 
breaking effects (and numerical effects). We sample the degeneracy surface at discrete points over a regular Q, u h 2 -h 
coordinate grid with parameter ranges chosen to give a good spread of \eS- 

In the absence of numerical effects, a Taylor expansion of C; in h and f2„/i 2 at fixed 9 should be accurate close to 
the fiducial model. The derivatives with respect to these parameters at the fiducial model can be already be inferred 
from the results in Sec. [IT] and Sec.|V] We thus expect physical differences in the spectra that are linear combinations 
of those described earlier. Note that our previous calculations have used either a flat code with massive neutrinos or 
a curved code with massless neutrinos. Our main aim here is to check the numerical robustness of the curved code 
with massive neutrinos. 

In Fig. [22j we plot the fractional differences between our degenerate models and the fiducial model with the former 
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FIG. 23: Contours of x^s over * ne two-dimensional degeneracy surface, 9 — const., in the Q, v h 2 -h-VlK space. We 
parametrize the surface by h and Q v h 2 and sample over a regular grid in these parameters. The contours enclose 
68% and 95% confidence regions. Left: default values of ISampleBoost, lAccuracyBoost and AccuracyBoost with 
high accuracy calculations (red dashed) and low accuracy (blue solid). Middle: accuracy parameters boosted to 2 
and high accuracy (red dashed) and low accuracy (blue solid). Right: high accuracy and default parameters (blue 
solid) and boosted to 2 (red dashed). The dashed (red) contours in the middle and right plots represent the true 
situation (i.e. with minimal numerical errors). Lensing and instrumental noise are not included. 



calculated at various accuracy settings. In all cases, the numerical errors are under control and for high accuracy with 
the other accuracy parameters at their default values the numerical errors are below the quoted 0.1% level. 

The physical differences revealed in the right-hand panels of Fig. [22] are consistent with our earlier findings. In- 
creasing the neutrino mass scale at fixed h, the changes in the curvature are rather small and the physical effects of 
the neutrino masses dominate the evolution seen in the power spectra. At fixed neutrino mass, the evolution in the 
spectra is due to the change in curvature, i.e. the late-ISW and projection effect discussed in Sec. ITT] 



In Fig. 23 we show contours of xts over t ne two-dimensional degeneracy surface of constant 9 parametrized by 
h and fi„/i z . The models are computed at various accuracy settings. The dashed (red) contours in the middle and 
right-hand plots have minimal numerical errors. The conclusions are similar to previous sections. Numerical noise 
in low accuracy calculations with default accuracy parameters distorts the likelihood significantly, shifting the best- 
fit point and artificially broadening the distribution in the f2„/i 2 direction. High accuracy calculations with default 
parameters perform better with only minor shifts in the contours. In this case, the high-frequency noise apparent in 
the middle column of Fig. [22] has little effect on parameter inferences. The physical degeneracy-breaking effects of 
curvature and neutrino masses have rather different spectral shapes (see Fig. 22 1 and this is reflected by the lack of 
correlation seen between parameters in Fig. |23| 



VII. MCMC SAMPLING ANALYSIS 



We have focused in this paper on maximizing the likelihood for different parameter values. This procedure has the 
advantage of being prior- independent, and is a rather stringent test of numerical accuracy because for convergence to a 
maximum of the likelihood the basin of attraction all has to be smooth and not give false numerical minima. However 
in practice most parameter analyses use sampling methods, which are rather less sensitive to small numerical errors. 
Figure 24 shows a set of MCMC samples from CosmoMC assuming a flat ACDM model, using the likelihood of 
Eq. Mm with Ci given by our fiducial ACDM model and Planck-like noise. The degeneracy that we explored in Sec. IV 



is clearly apparent, and by showing the value of other parameters for each posterior sample we can also clearly see 
the importance of other parameter variations in addition to fi m and H ; see Fig. |25| 

We used the high_accuracy_def ault setting for the MCMC runs, and obtained results that are fully consistent with 
the input parameters (as shown in Figs. 24 and 25). The numerical cost for high accuracy runs is quite manageable: 
assuming a good covariance matrix is available from previous forecasting runs, a ACDM six-parameter analysis using 
high accuracy settings only takes a few hours using two CPUs per chain (assuming the likelihood function is fast). 
By contrast, just increasing the accuracy parameters by hand to large values can easily waste very large amounts of 
computer time, and is unnecessary at Planck precision. 

It is worth noting that the parameter 9, which measures fairly directly the observed angular scale of the acoustic 
peaks [38] , is very well measured (at the 0.2% level); as such it is very sensitive to small changes in model, for example 
a change in the CMB temperature from 2.725 K to 2.726 K, or in the effective number of neutrinos from 3.04 to 3.046, 
will shift the calculated value of 9 by order la for fixed data. However the values of the other parameters are virtually 
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FIG. 24: MCMC samples from the posterior of Eq. ([6|) with Ci given by the fiducial ACDM model described in the 
text (marked with a cross) and Planck-like noise, including CMB lcnsing and polarization. The constraint on fl m h 3 
is nearly twenty times tighter than on the orthogonal direction. The degeneracy involves other parameters, so the 
accuracy, for example, of n s constraints from CMB alone is limited by this near degeneracy (see Fig. 25). 



unaffected compared to their error bars, since their relative precision is much lower because of the degeneracies. 

Small corrections due to changes in accuracy settings can easily be assessed from posterior samples by use of 
importance sampling, as explained in Ref. |33| . This only requires re-calculation of likelihoods at a subset of semi- 
independent samples from the original chains, and has the advantage of being almost independent of sampling noise 



(since the original samples and the importance-sampled samples are at the same points in parameter space). In Fig. 26 
we show fully-marginalized ID parameter constraints using various accuracy settings for the simplest six-parameter 
model; Fig. [27] shows equivalent results for non-flat and massive neutrino models. The default high accuracy setting 
produces significantly more accurate results than low accuracy, though even low accuracy only has a mean bias of 
about 20% of the error bar, and may be sufficient (at least as a base for importance sampling). The default high 
accuracy settings are still slightly biased compared to boosted accuracy results, though at a level that is < 10% of 
the error bar. The residual bias can partly be removed using lSampleBoost^ 1.2. 



A. More accurate multipole interpolation 



The main results of this paper have been calculated using the October 2011 version of CAMB, which interpolates 
the Ci in I using cubic-spline interpolation. The effect of the high accuracy setting (of Oct 2011) on this sampling is 
to increase the I sampling slightly (by ~ 20%). 

However the interpolation errors are in fact very well known, because we now have measurements of the CMB 
power spectrum that fix the location of the acoustic peaks very well. Differences in the C; between probable models 
only move the power spectrum by a small fractional amount compared to the best fit, so the interpolation errors 
evaluated for some fixed fiducial model near the best-fit model will also be very close to those in any nearby model. 
We can therefore simply subtract a fixed model of the interpolation errors; this makes the interpolation exact for the 
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FIG. 25: As Fig. 24 but exploring correlations between the other parameters. The parameter 9 is approximately 100 
times the ratio of the sound horizon to the angular diameter distance, and is very well constrained by the observed 
acoustic peak positions and nearly constant along the degeneracy [35]. Changes in fi m shift the angular scale of the 
CMB peak locations, though not very strongly, which is compensated by a shift in sound horizon at recombination 
due to a change in the physical matter density fl m h 2 [34] . However a change in physical matter densities changes 
the amplitudes of the acoustic peaks, so f^/i 2 and n s also have to change for partial compensation. 



fiducial model, very accurate for close models, and only increases the error for models that are very different (by an 
amount that goes to zero as the interpolation accuracy parameter is increased). This is equivalent to interpolating 
the difference between the Ci and a fiducial Cf , and adding the interpolated result to C{ calculated with dense I 
sampling 3 . This interpolation scheme has been implemented in the January 2012 version of CAMB, and numerical 
interpolation errors are compared in Fig. [28] The accuracy of the MCMC sampling results using this scheme is shown 
in Figs. 26 and 27 and is generally at a level < 5% of the error bar size with default high accuracy settings for Planck. 
This interpolation scheme is significantly more accurate for near-fiducial models, and as such high accuracy settings 
no longer need to increase the I sampling density; it is therefore faster than the high accuracy settings of October 
2011 used in the bulk of this paper (and higher accuracy; the sampling density could be further decreased slightly to 
gain additional speed advantage). 

We conclude that the default high accuracy settings should be adequate for analysis of the standard ACDM model 
at Planck sensitivity, and the January 2012 interpolation scheme achieves this at somewhat higher speed than using 
an interpolation scheme without error correction. Data analysis and physical model uncertainties are likely to be 
much larger than numerical issues. 



3 The interpolation is actually done on 1(1 + 1)Ci/2tt as this is more nearly constant on small scales. 
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FIG. 26: Left: Marginalized parameter constraints for the fiducial flat ACDM model with Planck-like noise, 
generated by using importance sampling from a well-converged MCMC run with default high accuracy settings (so 
all curves are generated using samples at the same parameter values). Solid black lines show the result using CAMB 

with default high accuracy settings; red dashed lines show the slightly biased result from using low accuracy 
settings. Magenta and blue curves use boosted high accuracy settings of 2 and 1.2 respectively, which are almost on 
top of each other. Right: Fractional difference in posterior means for different accuracy settings compared to the 
result from high accuracy with boosts set at 2 (as in the fiducial power spectrum). Default high accuracy settings 
are sufficient for biases to be < 10% of the error bar, but a small increase in ISampleBoost to ~ 1.2 is sufficient to 
remove most of the small residual bias (yellow bars show boosted accuracy) . Red bars show the result comparing to 
the new January 2012 version of CAMB at default high accuracy settings: this interpolates the difference between the 
Ci and a fiducial spectrum, significantly reducing interpolation errors for higher speed, so residual differences are 
then < 5% of the error bar. Less accuracy is likely to be required with more realistic data. 



VIII. CONCLUSIONS 

We have shown how parameter degeneracies remain very important even with high-precision CMB data. However, 
small degeneracy-breaking effects become quantitatively more important, especially CMB lensing on small-scales. The 
approximate acoustic-scale degeneracy in flat ACDM models will shrink significantly with future data, however it will 
remain a significant limitation on the ultimate precision of CMB-alone parameter inferences. 

We have shown that CAMB has sufficient numerical stability at high accuracy settings for degeneracies to be explored 
reliably, and quantified the size of residual numerical artefacts. Our tests conclude that CAMB is performing at or 
below its quoted numerical accuracy, and that the default high accuracy calculation is likely to be sufficient for 
interpretation of Planck data. Small residual biases due to multipole interpolation can be eliminated by subtracting 
the spline errors at a fiducial model, as implemented in the January 2012 version of CAMB; typical parameter biases 
are then < 5% of the random error expected with Planck. 
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FIG. 27: Effect of different accuracy settings on posterior parameter constraints for Planck-like noise compared to 

using high accuracy with accuracy parameters boosted to 2. Left: flat ACDM model with massive neutrinos 
(assumed degenerate); Centre: non-flat ACDM model with massless neutrinos; Right: non-flat ACDM model with 
massive neutrinos. As in the vanilla flat case, default high accuracy setting generally give results correct to within 
10% of the error bar, and using the new January 2012 interpolation scheme gives differences < 5% of the error bar. 
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Appendix A: Numerical evolution of massive neutrino perturbations 



The evolution of linear perturbations in the massive neutrinos is well understood |40j , but is more complicated and 
hence numerically slower than for massless species because neutrinos moving with different velocities evolve differently. 
Here we briefly describe CAMB's implementation as from July 2011. 

Scalar mode perturbations with comoving momentum q, multipole moment I, and wavenumber k evolve with 
conformal time as 



kv 
21 + 1 



[(I + l)0 l+ iF l+1 - LFf-i] + 



15 



Stati 



dlnF 
din (j 



0, 



where the time-dependent velocity is v = q/e and e is the comoving energy, and we follow the conventions of |T2"ll4"Tll4"2] . 
Here F without a subscript is the background distribution function, which is assumed to be Fermi-Dirac when the 
neutrinos are fully relativistic. For convenience we can divide through by defining vi = — 4jF;/(- dlrl ' 



kv 
21 + 1 



[Ivi- 



15 6 



When the neutrinos are relativistic so that v = 1, this is identical to the massless neutrino equation, and vi is then 
independent of q. The only difference with massive neutrinos is that they travel at a different (time dependent) speed 
once the mass becomes important. The hierarchies are truncated at Z ma x using 



u'i = kvvi-i — (I + 1) cotif vi. 
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FIG. 28: Power spectrum errors due only to cubic spline interpolation in I. The red and black lines on the 
right-hand side show the fractional error using a fixed I sampling with Al — 50 and Al = 42 respectively. The red 

line shows the error using the new January 2012 interpolation method based on interpolating (with sampling 
Al = 50) the difference to the fiducial model shown on the left (which is many sigma away from the model we are 
trying to interpolate). The new method is faster, and also significantly more accurate (ACi/Ci <C 10~ 3 ) even for 
models not very close to the assumed fiducial template; as the model gets closer to the template the interpolation 
error goes to zero. Similar results hold for the polarization, and using lensed/unlensed power spectra. 



Evolution is started when neutrinos are highly relativistic, when vi = J;, where J; are the massless neutrino 
multipoles. To get the leading correction from the mass we can write 

2 

171 A T 

n = Ji + ^aj ; 

2q z 

and use the series expansion v « 1 — a 2 m 2 / (2q 2 ) so that 

A J/ = ^ [iAJ,_! - p l+1 (l + l)AJ m ] - ^ [IJ,.! - p l+1 (l + 1)J, +1 ] . 

This lets us calculate the evolution of any q momentum mode while that mode has \am/q\ <C 1, which for light 
neutrinos and larger q can be a while: only one additional hierarchy has to be evolved (that for A J;) until modes 
start to become significantly non-relativistic, at which point the full mode equation must be integrated for a sample of 
momenta. The time saving from this approximation is not large in itself, but it does allow lower l max to be used when 
switching to integrating the momentum modes separately - without it the momentum modes need to be integrated 
from the beginning with roughly the same / max as the massless neutrinos. It also makes clear the leading q dependence 
of the perturbed distribution that we use below. 

CMBFAST, COSMICS and CAMB previously used a rather brute force integration sampling, using a fix grid of many 
q samples, which leads to numerically slow evolution: we would like to evolve as few different g-mode hierarchies as 
possible. Ref. [43 pointed out that one can do much better by using a more intelligent sparse sampling of q samples, 
e.g using Gauss-Laguerre quadrature. Here we describe a sampling scheme specifically optimized for the problem in 
hand. We need integrals of the form 

1 f°° g V w 
— I dq- -ttV vi 

in order to calculate the massive neutrino density, heat flux, and other perturbations enter the equations for the 
evolution of other species. In the perturbatively relativistic regime, so that we can do an expansion in a 2 m 2 /(2q 2 ) as 



2G 



before, the integrals are sums of terms involving integrals of the form 

[°° gV n 

At late times we also expect v ~ q, also giving terms roughly of this form, though the distribution has evolved away 
from anything simple. We do not attempt to integrate the distribution accurately in the intermediate sub-Hubblc 
regime where there can be oscillations in q: ignoring these seems to be harmless at required precision, presumably 
because averaged over time or k they are smoothed out. So the idea is to choose a sampling in q so that integrals 
with n = —4, —2. .2 are evaluated exactly, which gives a set of constraint equations for the points and weights that 
can be solved, and if more points are desired the solution can be made unique by adding other constraints or making 
a choice of a few points (e.g. at high q would expect Gauss-Laguerre point sampling to be nearly optimal 43] ). For 
3 points we find 4 the remarkably sparse sampling q = (0.913201,3.37517,7.79184) produces results accurate at the 
2 x 10~ 4 level with 

1 f°° a 4 e q v , 

lL ^(iW^?^ 

and kernel weights K = (0.0687359,3.31435,2.29911). A four-point sampling is accurate at the < 10~ 4 level, e.g. 
with q = (0.7, 2.62814, 5.90428, 12), K = (0.0200251, 1.84539, 3.52736, 0.289427). 

At late times further speedups are possible. For highly sub- horizon perturbation fcr>l (as for massless neutrinos), 
and once significantly non-relativistic, / max can be reduced down to 2 or 3 [IT] Once the neutrinos become very non- 
relativistic we can evolve velocity-integrated equations (i.e. a truncated fluid hierarchy); this is described in detail in 
Ref. [32] and was previously implemented in CAMB. 

Further speed ups may be possible using the fluid approximations of Ref. [23] , but the scheme described above is 
sufficient to dramatically reduce the computing compared to a naive evolution of a large sampling of q modes, so that 
calculations involving massive neutrinos only take 0(1) times longer than with massless neutrinos. 
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